Self-consistency and collective effects 
in semiclassical pairing theory 



V. I. Abrosimov a , D. M. Brink b , A. Dellafiore c '*, F. Matera c ' d 

11 Institute for Nuclear Research, 03028 Kiev, Ukraine 
O ■ 

CN ' b Oxford University, Oxford, U.K. 

\ c Istituto Nazionale di Fisica Nucleare, Sezione di Firenze 

Ph . d Dipartimento di Fisica, Universitd degli Studi di Firenze, via Sansone 1, I 50019 

£NJ \ Sesto F.no (Firenze), Italy 



Abstract 

c. 

A simple model, in which nuclei are represented as homogeneous spheres of symmet- 
ric nuclear matter, is used to study the effects of a self-consistent pairing interaction 
on the nuclear response. Effects due to the finite size of nuclei are suitably taken 
into account. The semiclassical equations of motion derived in a previous paper for 
the time-dependent Hartree-Fock-Bogoliubov problem are solved in an improved 
(linear) approximation in which the pairing field is allowed to oscillate and to be- 
come complex. The new solutions are in good agreement with the old ones and also 
with the result of well-known quantum approaches. The role of the Pauli principle 
in eliminating one possible set of solutions is also discussed. The pairing-field fluc- 
tuations have two main effects: they restore the particle-number symmetry which 
is broken in the constant-A approximation and introduce the possibility of col- 
lective eigenfrequencies of the system due to the pairing interaction. A numerical 
study with values of parameters appropriate for nuclei, shows an enhancement of 
the density-density strength function in the region of the low-energy giant octupole 
resonance, while no similar effect is present in the region of the high-energy octupole 
resonance and for the giant monopole and quadrupole resonances. 
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1 Introduction 



The study of pairing phenomena in nuclear physics has been an important 
topic for many years and different theoretical approaches have been developed. 
On one side, the seniority scheme, developed by Racah and collaborators [1,2], 
allows for the study of pairing in the shell-model framework. Modern devel- 
opments of this approach led to the so-called exact pairing formalism (see, for 
example, [3] and references therein). On the other side, the Bardeen-Cooper- 
Schrieffer (BCS) theory of superconductivity [4,5], originally developed for 
large, macroscopic systems, led Bohr, Mottelson and Pines [6] to point out 
the possible usefulness of the the BCS ideas for describing pairing phenomena 
in the finite nuclear systems. Further work by Belayev along this line showed 
that pairing correlations affect practically all low-energy properties of nuclei 
[7] . Since the early approaches, a great amount of work has been done on the 
effects of pairing in nuclei (see, e. g. [8] and references therein), however there 
are problems which deserve further attention, especially when finite-size ef- 
fects become important. One of these topics, which is of current interest both 
in the physics of nuclei [9] and of other mesoscopic systems [10], is the issue 
of possible collective phenomena associated with the pairing interaction. Here 
we address this problem by using a semiclassical approximation for the equa- 
tions of motion, but taking into due account the finite size of the many-body 
system. 

Recently, Urban and Schuck [11] have used a semiclassical approach to study 
the dynamics of trapped systems made of atomic fermions. These mesoscopic 
systems can contain a very large number of particles and can be viewed as 
laying between the macroscopic systems described by the BCS theory and 
the "small" nuclear systems. The semiclassical approximation of [11] leads to 
simplified equations, thus allowing for the solution of problems that would 
become unmanageable in a fully quantum approach for a large number of 
constituents. 

The semiclassical approximation is valid for small or mesoscopic systems with 
a characteristic size R, provided ppR/h » 1, where pp is the Fermi mo- 
mentum. This is equivalent to the condition that hu « e F , where Uq is a 
frequency characteristic of the classical orbits and ep the Fermi energy. We 
consider systems which, in the absence of pairing, are described by the hamil- 
tonian h = p 2 /2m + Vo(r). Except for special cases, the classical orbits in the 
potential Vo(r) are not exactly periodic, they are rather multiply periodic (see, 
for example, [12] p. 463), leading to eigenfrequencies u n = J2 a n a^a , which 
are harmonics of the fundamental frequencies oo a ~ uq [13]. In nuclei, these 
eigenfrequencies are related to the frequencies of giant resonances . 

Since in heavy nuclei the condition p F R/h » 1 is well satisfied, in [14] we 
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have used an approach similar to that of Ref. [11] to study the linear response 
of (heavy) nuclei, with the aim of developing a simplified tool for the study of 
pairing effects in low-energy nuclear excitations. In that paper a drastic ap- 
proximation has been made: the pairing field A(r, p, t) has been approximated 
with a constant phenomenological parameter A (both in the statics and dy- 
namics). It is well known that such an approximation violates the continuity 
equation and consequently introduces spurious contributions into the density 
strength function and its energy- weighted sum rule (EWSR). This problem 
was solved in [14] by means of an appropriate prescription which enforces 
upon the density fluctuations the constraint coming from the continuity equa- 
tion, however a prescription is not very satisfactory from a theoretical point 
of view and here we want to improve on that approach by introducing better 
kinetic equations that do not need any prescription in order to satisfy the 
continuity equation and the EWSR. The resulting equations allow also for the 
study of collective effects associated with the pairing interaction. 

Both here and in Ref. [14], pairing is included in semiclassical equations of 
motion as in Ref. [11]. We focus on weak pairing where the pairing gap A << 
hu>o- This paper is concerned with the frequency response of nuclei and in our 
formalism the pairing gap is a frequency. In the limit of small h the gap may 
be small, while the frequency A/h remains finite. Mesoscopic systems with 
A >> hujQ can also be studied with the formalism of this paper, but they are 
not discussed here. 

While in [14] the (complex) time-dependent Wigner-transformed pairing field 
A(r, p, t) was approximated with a real, constant, phenomenological parame- 
ter A, here we study the improved approximation in which, for small density 
fluctuations 



A(r, p, t) = A (r, p) + <JA(r, p, t) « A + 5A r (r, t) + i5A\y, t) . (1) 



Thus, in the present approximation, the pairing field is allowed to oscillate and 
to become complex, the possible momentum dependence of the complex fluc- 
tuations is neglected though, in order to simplify the theory. The static pairing 
field A (r, p) is approximated with the same phenomenological constant used 
in [14], hence our approach is not fully self-consistent, however, the pairing- 
field fluctuations SA r (r,t) and 5A l (r,t) are derived from self-consistent rela- 
tions. We study the new set of equations of motion that arise in this improved 
approximation and determine the effects of the additional terms on the density 
response function of nuclei. 
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2 Basic equations and approximations 



2. 1 Equilibrium properties 



We assume that our system is saturated both in spin and isospin space, so 
we do not need to introduce explicitly the spin and isospin variables and 
the physics can be discussed, semiclassically, in terms of the two equilibrium 
phase-space distributions p ( r iP)) an d K (r,p) which, according to Ref. [15], 
p. 550, are given by 1 



lA (r,p) 
with the quasiparticle energy 



E(r, p) = JA%(r, p) + (h (r, p) - . (4) 



The chemical potential \i is determined by the number of nucleons A through 
the relation 

4 r 

A = To tvs / drdpp (r, p) . (5) 



(2nhf 



The equilibrium hamiltonian 

Mr, P ) = f^ + VUr) (6) 



contains the (Hartree) mean field Vo(r), which should be evaluated self-consistently, 
however in the following we use a phenomenological potential well instead. Like 
in [14], we approximate the static nuclear mean field with a spherical square- 
well potential of radius R, this choice allows us to take into account finite-size 
effects and, at the same time, to recover the simplicity of homogeneous sys- 
tems. 



In the following, we also approximate the equilibrium pairing field A (r, p) 
with the phenomenological parameter A, which, in heavy nuclei takes values 

1 note the opposite sign of our function k , compared to the function n in Refs. 
[11] and [15] 
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between 1 and 1.5 MeV ([16], p. 170). As a consequence of this approxima- 
tion, both the static and dynamic equations are considerably simplified. The 
static distributions, in particular, become a function of the particle energy 
e = h (r, p) alone, like in an infinite homogeneous system: 



Po(e)4(l-^), (7) 



E(e) = J A* + (e - , (9) 

Po[e) ~ de ~ 2ES( C )' (iUj 
_ dnoje) _ A(e-/i) 



2.2 Dynamical equations for small amplitudes 
The dynamical equations of motion are [11,14] 



ik^- = ih{h, p} - 2iIm[A*(r, p, t)n] + ^Re{A*(r, p, *),«}, (12) 
ih— =2{h-n)K- A(r, p, t)(2p ev - 1) + ih{A(r, p, t),p od } , (13) 

with 



Pev = ^[p{r,p,t) + p(r,-p,i)] (14) 
Pod = ^[p(r,p,t) -p(r, -p,f)] • (15) 

The equations of motion are gauge- invariant in a self-consistent theory, be- 
cause multiplying both k(t, p,t) and A(r, p, t) by a common phase factor 
exp(ix) leaves the equations unchanged. Here we make a particular choice of 
gauge where the equilibrium quantities n and A are both real. In the small- 
amplitude limit, we write the equations of motion for p and k in terms of their 
deviations from equilibrium. To first order, 5k t gives the change in magnitude 
of k and 5k 1 k$x is proportional to the change x m the phase of k. 

Thus, we assume that, at time t — 0, the system is perturbed by a weak 
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external driving field of the kind 
5h(r,t) = r)6(t)Q(v), 



(16) 



where 5(t) is a Dirac 5-function in time and n is a parameter specifying the 
strength of the external field, then at time t > 0, we have to deal with the 



quantities 

p-(r,p,t) = p ( e ) + ^(r,p,t), (17) 

p od (r,p,t) = 5p od (v,p,t), (18) 

k(t, p, t) = K (e) + 5K r {r, p, t) + i5n % (r, p, t) , (19) 

A(r, t) = A + 5A r (r, t) + iSA^r, t) . (20) 

In the present improved approximation, the dynamic equations derived in [14] 
become (cf. Eqs. (43-46) of [14]) 

ihd t 5p ev = ih{h , 5p od } - 2iA[SK i + -^-] (21) 

ihd t 5p od = ih{h , 5p ev } + ih{5h, p }+ih{5A r , k } , (22) 

5A r 

-hd t 5K l = 2(e - aO [5k t + — — ] + 2K 5h - 2A5p ev , (23) 

2E(e) 

5A { 

hd t Sn r = 2(e - p) [Sk 1 + ——] . (24) 



In Ref. [14] the fourth equation of motion (EOM) has been replaced by a 
supplementary condition enforced by the Pauli principle (cf. Eq. (54) of [14]). 
This was done in order to simplify the formalism. Here instead, we use the 
four EOM (21-24) as our staring point and will show that the two methods 
give similar results. 

In Eqs. (21-24) terms of order h 2 or higher, have been neglected. 
2.3 Self- consistency and continuity equation 

Compared to [14], the present approach contains the two extra unknown func- 
tions 5A r (r, t) and 5A l (r, t), hence we need two additional equations in order 
to determine these quantities. 

In a self-consistent theory, the changes in the pairing field are related to the 
changes in the anomalous density. There are many possible choices for the self- 
consistency relation. Two minimum requirements are that the total particle 
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number should be conserved and that the particle density and current density 
should satisfy a continuity equation. A choice that satisfies both requirements 
is obtained from the self-consistency relation (gap equation) written in the 
form (cf. Eq. (11) of [11], see also [18]): 



7*( 



A(r,t) p 2 /m 



) = 1. (25) 



Here g is a parameter determining the strength of the pairing interaction. 
We have assumed that the p-dependence of the dynamic pairing field can be 
neglected. By differentiating Eq. (25), we get the first-order relation 

J dp{5K(r,p,t) - K (r,p) 5 -^^-) =0, (26) 

where K (r, p) and A are real equilibrium quantities, while 6k and 5A are 
their complex fluctuations. 

The real and imaginary parts of Eq. (26) give the two independent relations: 

J dp(y (r, p, t) - «o(r, P)^£^) = , (27) 
J dp(5K i (v, p, t) - «o(r, P)^^) = . (28) 

These conditions, based on Eq. (25), take into account the residual pairing 
interaction in a self-consistent way. Note that, since the conditions (27, 28) 
do not depend on g, they are valid both for weak and strong pairing. For our 
purpose, we do not need further information about the pairing interaction. 

Using Eq. (8) gives 



Jdp(8K'(r,p,t) + ^^)=0, (29) 
/rfp(^(r,p,t) + ^^)=0. (30) 

Self-consistency conditions similar to these have been used by the authors of 
Ref. [10] in their quantum calculations for infinite homogeneous systems. 

It is important to check that, even with the approximate form (29, 30) of the 
self-consistency relations, the continuity equation is still satisfied. 
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Integrating Eq. (21) over p gives 



ih dtj dp5p ev (r,p,t) = thj dp{h ,5p od } - 2tA J dptfj + ^y] .(31) 



2E(eY 



The last integral vanishes because of Eq. (30), giving 

d t 5g(r,t) = -d r ■ j(r,t) + J dpd r h ■ d p 5p od , (32) 

with the density fluctuation 

Mr, t) = J dpSp ev (v : p, t) (33) 

and the current density 

j ' r ' t »=(2^/*^< r - P ' t »' < 34 » 

For a Hamiltonian of the kind (6), the integral in Eq. (32) can be shown to 
vanish and the continuity equation is satisfied. This is a crucial difference 
between the present approach and that of [14]. 

3 Fourier expansions 

First of all, we take the Fourier transform in time of the various time-dependent 
fluctuations in Eqs. (17-20). Since these quantities are nonvanishing only for 
t > 0, we suppose that u has a vanishingly small positive imaginary part ie, 
to insure convergence of the integrals f(u) = $™ 00 dte v * t f(t). In the following, 
the imaginary part of u will not be written explicitly. 

Then, following [14], we introduce the following Fourier expansions based on 
the method of action- angle variables (<fr are the angle variables, while I are 
the action variables): 

5h{r,u) = Y,^h n {I,u)e in ^ , (35) 

n 

5p(r,±P,^) = E^n(I^)e m -* ) (36) 

n 

Sp e n v = \[Spt + S Pn ] (37) 
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&tf = \[&pt-&PZ\ (38) 

< 5« r ' i (r,p,o;) = X;*«5 i (I,w)e in -*, (39) 

n 

5A r >\r, w) = £ 5A^(I, w)e m '* . (40) 

n 

Moreover we use the relation 

{/,M = E*«W» e<n '*» ( 41 ) 



where a; n = ^ Q n a w a are the eigenfrequencies of the uncorrelated system 
[13], to obtain from the dynamic equations (21-24) the system of algebraic 
equations 



hu5p e : = hw n 5p° n d - 2iA[SKl 1 + , (42) 

/lo;5p° d = hcu n 5pn - hu n p' 5h n - hu n K' SA r n , (43) 

ihwSKi = 2(e - /i) + ^] + 2k ^„ - 2 , (44) 

-^< = 2( e - / i)[5< + ^]. (45) 



From Eqs. (42, 43) we get 2 

( W 2 - col) 5p e : = -u 2 n p' 5h n - co 2 n K' 5K - 2^A5< - , (46) 



while Eqs. (44, 45) give 

, i u» , 5A l n . 5A r n . 2A_ ev . 2k 0£1 , 
S < = 2 a K-FT _ """FT + w W - <ta„ . 47 



By combining the last two equations, we obtain 



5A i 5A r 5h r 

5K,' t n = A n (u)—^ - iuu^AJu)-^ + 2iuAA n (u) ' 



(48) 



with 

2 From now on, we use units h = c = 1. 
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2 n2 a A 2 W n 



D n (w) = w 2 - - 4A 2 ^-^ = W 2 - n 2 - 4A 2 ^p^ , (50) 
and 

il = 2E = ^4A 2 + w 2 , (51) 

^ = 2(e-/i), (52) 

^ = ^ + Q 2 . (53) 

In a similar way, we get also 

S< = ^[1 + + *4nM^ - 2Au V A n (u;)^ (54) 



and 



^2 _ 2 



By using the last two equations, it can be easily checked that, for the mode 
n = (0, 0, 0), the following relation is satisfied: 

*«» = • (56) 



In [14] the equation of motion (45) was replaced by this relation, which is 
a condition required by the Pauli principle. Thus we see that, for the mode 
n = 0, the Pauli principle constraint (56) is authomatically satisfied by the 
solution of the equations of motion (42-45). The more general case n ^ is 
not so immediate and it is discussed in the Appendix. 
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4 Constant-A approximation 



The first term in Eq. (55) gives the constant-A part of the density response. 
By using Eq. (53), it can be written as 

, uiiu 1 V 1 ) . 4A2 . (u 2 V 1 ) \ 4A 



The first part of this expression contains a factor l/Q 3 which, as a function of 
e, has a very narrow peak at the Fermi surface when A is small compared to 
li. The second part instead contains a factor l/Q, which is also peaked at the 
Fermi surface, but is more spread out. This second term gives an unwanted 
contribution to the energy- weighted sum rule [14]. In the following we will 
check explicitly that the self-consistent pairing field fluctuations <5A„ restore 
the sum rule to its correct value. 

In the constant-A approximation the eigenfrequencies of the system are de- 
termined by the poles in Eq. (57), that is, by the condition 

(u 2 - u 2 n )D n (u) =u A - oolu 2 + u>lu>l = , (58) 
with solutions 

: (n,I) = k(I)(l± 



(59) 



These expressions look rather unusual, however they are closely related to the 
eigenfrequencies found in [14] and also, as will be shown in the following, to 
the energy of two-quasiparticle excitations. 

First of all we notice that, when n = 0, the eigenfrequencies u + coincide with 
the eigenfrequencies u n found in [14] since 

c^(n = 0,I) = fi 2 (e), (60) 



while 

oo 2 _(n = 0,1) = 0. (61) 
Moreover, for the modes n ^ 0, we have that 
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"l(n,I) 
co 2 _(n,I) 



= ^(I), 



0. 



(62) 
(63) 



As a general rule, the range of particle energies which is most important is 
that within a distance of a few times A from the Fermi surface /x, hence is 
of order A and the dimensionless parameter 



a n (I) 



(64) 



both for small (w n (I) >> A) and large (w n (I) << A) systems. In these cases, 
by expanding the square root in Eq. (59), we have that 



(n,I) 
(n,I) 



^(1) 



(65) 
(66) 



This is an excellent approximation in the first term of (57) because of the factor 

2 2 

^j^y- It is interesting to note that the factor (u 2 — < ^ E ) in the numerator of 

5p^ approximately cancels the poles (u 2 — uj 2 _ (n, I)). The same cancellation 
occurs in the coefficient of 5A^ in Eq. (55), while it does not occur for the 
coefficient of SA % n . 

In conclusion, in the constant- A approximation, the solutions of the equations 
of motion (42-45) have two branches, corresponding to the eigenfrequencies 
w+(n, I) and w_(n, I). These solutions do not satify the constraint (56) exactly 
when n^O, but, as shown in the Appendix, one branch does, approximately. 
Since a n (I) << 1, the eigenfrequencies of the two branches are approximated 
by uj\ ~ w„ and uj 2 _ ~ u^u 2 /^. The poles corresponding to w_ are approxi- 
mately canceled by the term (u 2 — u^u^/u^m the numerators of Eq. (57). 

Thus, the constant-A part of the present approach, which is based on the four 
equations of motion (42-45), gives results in good agreement with those of Ref. 
[14] in which the fourth equation of motion was replaced by the constraint (56). 



4-1 Relation with two-quasiparticle excitations 



The eigenfrequencies (59) are related also to the energy of two-quasiparticle 
excitations. 
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In the quantum BCS theory , the simplest excitations of an even nucleus are 
two-quasiparticle states, with energy (cf., for example, Sect. 3.3 of [8]) 



S(e p , e h ) = E P + E h = yj(e p - /i) 2 + A 2 + - e h f + A 2 . 



(67) 



If we change variables from (e p , e h ) to (e p h,e), with 




(68) 
(69) 



then 



3 = J + (c - »)Y + A2 + V " ( e ~ " ^ + A2 • 



From the work on normal systems [13], we know that the eigenfrequencies u n 
give a semiclassical approximation to the particle-hole energies 3 



The last relation is exact, as can be easily checked by squaring both sides. 

Thus, the eigenfrequencies w+(n, e) correspond to two-quasiparticle excita- 
tions. There are two approximations involved: one is the relation (70), and 
the other is the replacement e — > e in the particle energy. 

In the same way it is possible to show that the eigenfrequencies w_(n, e) corre- 
spond to the combination E p — E^. Anderson [21] called 'unphysical' the solu- 
tions corresponding to this combination of quasiparticle energies and Valatin 
[22] pointed out that these solutions are eliminated by the supplementary 
condition required by the Pauli principle (see Appendix). Here we find that 
the poles corresponding to these 'unphysical' eigenfrequencies are practically 
canceled by corresponding zeros in the numerator of the density fluctuations 



3 In order to simplify the notation, we will write also oj n (e), instead of w n (I) 



tph ~ w n (e) , 



(70) 



then 




(71) 



(57). 



c 
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5 The self-consistency relations 



The initial system of coupled differential equations has been reduced to a 
linear algebraic system, which is supplemented by the two integral relations 

/jpfc(.,M)+ a ^;" ) )= l - < 72 > 

fdp(S^T,p, U ) + S -^f)=0, (73) 

which are the Fourier transform in time of Eqs. (29, 30). 

Multiplying both sides of these relations by e - " 1 '* and integrating over r, we 
get 



6A i (r,u)\ _ in .^ 



dr J dp(6K*(r,p,u)+ 2 M e ) ) e =°> ( 74 ) 
Jdxj dp(«W(r, p, u) + ^Effi *'* 1 '* = ° " (?5) 



Now, changing variables from (r, p) to (I, and using the orthogonality of 
the functions e m $ , we have 



/fl(K(I, w ) + ^g^) = 0, (76) 
/dl(K(I,o;) + ^|^)=0. (77) 



Finally, by using the four equations (42-45), the two integral relations (76, 77) 
can be written in the form 



J dI[a 11 (n,I,a;)<JA;(I,a;) + a 12 (n,I,a;)(JAi 1 (I,a;)] = J c?I&i(n, l,uj)5h n (T) , 
| dl[a 2 i(n, I, w)5A^(I, w) + a 2 2(n, I, w)5A* n (I, u)] = J dlb 2 (n, I, u)5h n (I) , 

with 

w„ ujn / 4A 2 w 2 \ 

a n (n,I,a;) = -^-A n H = -^^-^ (l + T^TJ^j > ( 80 ) 
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1 + A n (w) u 2 



oo 2 1 

a 21 (n,I, W ) = -|A n H + - (82) 



fi 3 nv 1 n 



w 2 -4A 2 (l + 



/ j \ c^l + A n (w) . 



7 / t \ 2icuAA n (w) 2iwA , 4A 2 w 2 x 



fi 3 "~ OD n (u;)V VI 2 oo 2 -ul 

2AuyA n (a;) 2Aay , 4A 2 a; 2 \ 

62(n ' w) = — ns — = twjMjv- + T^^z) • (85) 



(84) 



The external field fluctuations (16) depend on position and time, but not on 
the momentum p. For a uniform system the self-consistency relations lead to 
pairing- field fluctuations 5A(r,u) which are also momentum-independent and 
Eqs. (78, 79) become separable, giving an algebraic system. For finite systems, 
Eqs. (78, 79) are not separable because the coefficients Sh n (l) can depend on 
the action variables I. 



6 'Ansatz' and solution 



The pairing-field fluctuations are given by the solution of the two coupled 
integral equations (78, 79). In the spirit of a simplified approach, here we 
make an Ansatz that simplifies the problem: we assume that 

5A r n (I,u) = F:(Lu)5h n (I), (86) 
5A i n (I,u) = F^u)5h n (I), (87) 

with the functions F^ l (u) given by the solution of the algebraic system 



F n » J dIa n (n,I,u;) + F>) f dla 12 (n, — J dlh(n,l,uj) , 
F n » J dIa 21 (n,I,u) + F^uj) J dla 22 (n, I, u) = J dlb 2 (n, I, u) . (89) 

Then, defining 

A ij (n,u) = J dlaij(n,l,uj) , (90) 
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B 

we have 



(n,w) = J dI6i(n, I,a;)) , 



(91) 



B\A 22 — B 2 Ai2 

Ai\A 2 2 — A 2 \A\2 ' 
B2A11 — B1A21 

A11A22 — A21A12 



(92) 
(93) 



The six integrals Ay and Bi are conveniently expressed in terms of the four 
basic integrals 



as 



I 2 (n,u) = J 



Q(e) D n (oo) 
1 1 



1 w 2 



« 3 (e) D n (u) oji-ui 



A n = -iu(h + 4A 2 / 3 ) = -iul[ , 



A 12 = oj 2 I 2 , 



A 2l = {u 2 - 4A 2 )/ 2 - 16A 4 / 4 = oo 2 h - 4A 2 I 2 
A 22 = iuh , 

Bx = -iu2A(I 2 + 4A 2 / 4 ) = -iu2AI' 2 , 
B 2 = 2A(J 1 + 4A 2 J 3 ) = 2A/(, 



(94) 
(95) 
(96) 
(97) 



(98) 
(99) 
(100) 
(101) 
(102) 
(103) 



with 



J( = / 1 + 4A 2 / 3 , 
7 2 = / 2 + 4A 2 / 4 . 

In terms of these integrals, the solution (92, 93) reads 



(104) 
(105) 



F n » = 2A 

2A r 

F» = - 1 + 4A 

VjJ L 



4A 2 (/ 1 / 4 - I 2 I 3 



hI[-I 2 (u 2 I 2 -4A 2 I 2 



I 3 I[ - h{u 2 I 2 - 4A 2 / 2 



J 1 /(-/ 2 (a; 2 / 2 -4A2^) 



(106) 
(107) 
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With these expressions, the density fluctuations and the density-density re- 
sponse function can be evaluated explicitly. This is what we do in the follow- 
ing, by neglecting that part of the mean-field fluctuations that is not related 
to pairing, hence the expressions derived in the following correspond to the 
static mean-field approximation, called 'zero-order' in Ref. [13]. Within this 
approximation, the fluctuations of the hamiltonian are 

5h n {I, u>) = r]Q n (I) = f d$e- in *Q(r) . (108) 



6. 1 The mode n = 



The mode n = (0, 0, 0) is particularly interesting because, as shown in [14], this 
mode is the only one giving a spurious fluctuation of the number of particles 
in the system. When n = 0, the two integrals I 3 and J 4 vanish, I[ = I ly I' 2 = I 2 
and the solution is expressed in terms of the integrals I 1}2 alone, giving 



^n=oH = 0, (109) 
2A 

n =0 M = — , (no) 

LUJ 

with the possible exception of points where the denominators in Eqs. (106,107) 
might happen to vanish. 

By using the result (109, 110), it can be easily shown that the spurious con- 
tribution to the mode n = in the density fluctuations of the constant-A 
approximation, is exactly canceled by the fluctuations of the imaginary pair- 
ing field, since Eq. (55) gives 



fl 2 4A 2 1 4A 2 



Another interesting feature of the mode n = is that the condition 

A u (n, lu)A 22 (ii, u) - A 21 (n, u)A 12 (n, u) = , (112) 



that determines the possible collective eigenfrequencies of the system, acquires 
a particularly simple form and does not depend on the equilibrium mean field: 

^ViM-(^-4A 2 )/f(u;)] = 0, (113) 



with 
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'^i/^m^m' (114) 

U^J^Wj rf-Wie) - (U5) 

By evaluating these integrals, we find that the dispersion relation (113) has 
only the solution u 2 = 0. This zero-frequency solution corresponds to the 
Anderson-Goldstone-Nambu mode, associated with rotations in gauge space 
(see e.g. Ch. 4 of [8]). 

The other possibility 

[J>)-(w 2 -4A 2 )if(w)] = (116) 



instead, is reminiscent of the dispersion relation of pairing vibrations (cf., e.g. 
Eq. (5.23) of [8], note the correspondence between our integrals hpi^) and 
the quantities B and A of Ref. [8]), however we find that the relation (116) 
has no solution in the interval u = — 30 MeV. The quantum counterpart of 
Eq. (116) instead, has many solutions, corresponding to single-particle levels 
6j close to the Fermi surface (cf. e.g. Eqs. (J. 31) and (J. 32) of [8]), but these 
solutions have a single-particle, rather than a collective character. 



7 Response function and dispersion relation 



Equations (55, 92, 93) allow us to calculate the density response function, 
defined as 

K(uj) = - [ drQ*(r)Sg(r,uj), (117) 



where 5g(r, ou) is the time Fourier transform of the density fluctuation (33). 
According to Eq. (55), the density response function is given by the sum of 
three terms: 

k(u) = n cd {u) + n r {u) + n i {u), (us) 



with the first term containing the contribution of the constant-A approxima- 
tion, the second and third terms, the contributions of the real and imaginary 
parts of the pairing-filed fluctuations, respectively. The three components of 
the density response function are given by 
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ft c V) = 4E / dIC?(I,u)Q n (L)Q* n (L), (119) 

n •* 

n r {u) = ±Y, F » I d\Clil,uj)Q n {l)Qlil) , (120) 

n J 

7^H = 4£F>) / dIC;(I,a;)Q n (I)g;(I), (121) 



with 



x w^-w^ 4A 2 

c » (I ^»=(^I)t(t)^ (122) 

w 2 w^ — w 2 w 2 2Aw„ 

^"'' l^wi V' (123) 

^ ") = -«-' < 124 > 

The solutions (106, 107) can be written as 



F n » = 2AF», (125) 
2A 2A ~ • 

^H = — + — KM, (126) 



with 



nl J 1 /(-/ 2 (o; 2 / 2 -4A^)' 



(127) 
(128) 

Since the function F^(u) is the sum of the two terms in Eq. (126), the response 
function (121) also becomes the sum of two response functions 



2A r 

k\uj)=aY j — / rfic;(i,u;)Q n (i)g;(i) + ^H, (130) 

n iuj J 

with 

^H= 4 E^nM / cnci(i, w )Q n (i)o;(i). (131) 



19 



The first part of TZ l (u) cancels the spurious part of lZ cd (u) (the term ou 2 Q 2 of 
the product w 2 w 2 = w 2 w 2 + oj 2 VL 2 in the numerator of (122)) and the density 
response function becomes 

K{u) = n cd (co) + n r (co) + n'iw) , (132) 



This response function can be written in a more compact form as 



ftM=4£ / dlG n (I,u)Q n (T)Q* n (I), 



(133) 



with 



AA 2 ( u^-ug 



( W * - u£)£> n (a;) 



(134) 



w 2 ft 2 



The response function (133) can be considered as a finite-system version of 
the response function derived by the authors of Ref. [10] for uniform systems 
(cf. Eq. (B22) of [10]). The dispersion relation of the Bogoliubov-Anderson 
mode of [10] is here replaced by the dispersion relation 



Ji(n,u;)Jj(n,u;) - J 2 (n, w)[w 2 / 2 (n, w) -4A 2 ^(n,w)] = (135) 

that takes into account the finite size of the system. This dispersion relation 
is the main result of our paper and in the next section we will study the 
possibility of a collective mode, analogous to that studied in [10], in 'small' 
systems like nuclei. 

It can also be easily checked, simply by looking at the large-cu behavior of G n 
[20], that the EWSR is satisfied by our response function (132), because the 
spurious contribution of 1Z cd (u) is canceled by the pairing-field fluctuations 
(more specifically, by the first term on the right-hand side of Eq. (126)). 



8 Spherical cavity 



The four three-dimensional integrals (94-97), which appear in the dispersion 
relation (135), are not simple when n ^ 0. In a spherical systems, however, 
these integrals can be reduced to two-dimensional. In such a system, the action 
variables are conveniently defined as ([12], p. 476, see also [19]) 
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■Z ) 



drp. 



(136) 
(137) 

(138) 



where A is the magnitude of the particle angular momentum, X z its ^-component 
and p r the radial component of the particle momentum. It is convenient to 
change variables from I to (e, A, cos/3), with cos/3 = A 2 /A, the Jacobian of the 
transformation is 



where T(e, A) is the period of radial motion. Because of the spherical symmetry, 
the eigenfrequencies w n (I) do not depend on X z . Moreover, as shown in [13], 
the eigenfrequencies can be labeled with only two components of the integer 
vector n. We put n 3 = n, n 2 = N and n\ becomes redundant. Thus, in 
spherical systems 

w n (I) ->• u nN (e, A) « e r +n,i+N - £ r ,i (140) 

(in a spherical mean field, the single-particle levels e r ; depend only on two 
quantum numbers). 

Also, by expanding the external field Q(r) in multipoles as 



J(e,A) 



AT(e, A) 
2^ 



(139) 



Q{v) = Y J Q L {r)Y LM {v), 



(141) 



LM 




(142) 



L 



with 




x|gS( e ,A)| 2 G(iV,n,6,A,a;), 



(143) 



and 




(144) 
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^F r (N,n,u) 



W 4 - U 2 U 2 lN + ulul N 



F l {N,n,u 



» 



u 4 - u 2 u 2 nN + u 2 u 2 nN 



In Eq. (143), only terms in which the integer iV has the same parity as L appear 
in the sum over N, otherwise the spherical harmonics Ylat(|,0) vanish. The 
Fourier coefficients Q n (l) in Eq. (133) are replaced by the radial coefficients 
<5^jv(e, A), defined as in [13], which correspond to the radial matrix elements 
of the quantum approach. 

The quantity A(e) is the maximum possible value of A for a particle with 
energy e. 

For a spherical cavity of radius R, A(e) = \j2meR and the eigenfrequencies 
cu n 7v(e, A) are conveniently written as 

w„Ar(e, A) = u F {e)s nN {x) , (145) 

where x is a dimensionless parameter, which is related to the particle angular 
momentum by 

x = sin a, (146) 



A 



(147) 



cos a = 



A(e)" 



The functions s n ]y(x) are 



SnN(x) 



nit + N arcsin(x) 



(148) 



x 



while the frequency up(e) is given by 





(149) 



R 



The multipole response function (143) becomes 
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x j dxx 2 \Q { ^(x)\ 2 G(N,n,e,x,uj). 
o 

The functions G{N,n,x,oS) are given by Eq.(144), with 

h(N,n,u) 
I 2 (N,n,co) 
h(N,n,u) 
h{N, n, u) 
Ji(N,n,e,u)) 
J 2 (N,n,e,u) 

and Af is a normalization factor which cancels in the final expression. 
The denominator D(N,n,e,x,u>) is given explicitly by 

Dl N, n, e, A, tf( t ) - f**^ ■ 

The dispersion relation (135) takes the form 
K(N,n,u) =0, 

with 

if (iV, n, w) = ii(7V, n, w)i((7V, n, u) 

-I 2 (N,n,u)[uj 2 I 2 (N,n,uj) - AA 2 I' 2 (N,n,u)} . 
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oo 



= W y* deye^-yJi(A^,n,e,u;) , 
= Afjde v Te ^u 2 F (e)J 2 (N, n, e, w) , 

oo ^ 

= M de^——uj 2 F (e)J 2 (N,n,e,uj), 
J fi 3 (e) 



o 

dxX 2 —— r , 

o D(N,n,e,x,u) 



i 

dxx 



o 



1 1 S nN( X ) 

D(N,n,e,x,uj)uj 2 - uj 2 F (e)s 2 nN {x) ' 



9 Results 



The first figure shows the absolute value of K(N, n, id) for a few low-energy 
modes (N, n) in the range id = — 30 MeV for a spherical cavity containing 
A = 208 nucleons within a radius R — 1.2 A 1 / 3 fm. The chemical potential \i is 
determined from Eq. (5) to be fi ~ 33.33 MeV and a value A = 1 Mev has been 
assumed for the static pairing gap. Except for the modes (N = ±l,n = 0) 
(top right panel), all other modes display a smooth behavior suggesting that 
for these modes the dispersion relation (158) has no solution in this range of 
u. The modes N = ±l,n = instead, have a minimum at u ~ 7.7 MeV, im- 
plying that, for this value of id, the dispersion relation (158) is approximately 
satisfied. Thus we can expect (mild) collective effects in the odd multipole 
channels, to which these modes contribute. The dipole response gets its main 
contribution from the mode (N — ±1, n — 0), however this channel is spu- 
rious because of the breaking of translation invariance, so the first physically 
meaningful collective effect is expected in the low-energy octupole resonance 
(corresponding to L — 3, N — ±1, n — 0) while the high-energy octupole res- 
onance (corresponding to L = 3, N = ±3,n = 0) is not affected by collective 
effects (cf. bottom right panel of Fig.l) 



Fig. 1. Absolute value of K(N,n,u) as a function of lo (in MeV) for n = and 
N = (top left), N = ±1 (top right), N = ±2 (bottom left), N = ±3 (bottom 
right). The curve corresponding to the modes N = ±1, n = has a local minimum 
at co « 7.7 MeV, corresponding to an approximate solution of the dispersion relation 
K(N,n,u) = 0. 

The second figure shows an enlargement of the top right panel of Fig. 1, in 
the region of the minimum. 
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Fig. 2. Absolute value of the function K(N = 1, n = 0, oj) in the region of the local 
minimum. 

The other figures (Figs. 3-6) show the strength function 



S(u) = lmK L (uj) 

7T 



(160) 



associated with the multipolarities L = 0, 2, 3. 

The radial Fourier coefficients needed for the monopole and quadrupole re- 
sponse functions are (for the quadrupole and octupole response, we take 
Ql(j') = f' L ■, while for the monopole response, we take Ql(j') = r L+2 ) 



Q^(x) = (-) n R 



n d2 



S no( X ) 



(n^0) 



Q ( nUx) = {-TR 2 - 



[1 ±2 



X) 



(161) 

(162) 
(163) 



while those entering the octupole response are 



X 



-) n R 3 -r 



7 nN 



[1 + -N 



+ 4(|7V| - 1) 



3" s nN (x) 
l-x\ 



S un( x ) 



N = ±1,±3. 



(164) 
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Fig. 3. Contribution of the modes (N = 0, n = ±1) to the monopole response 
function of a spherical cavity containing A = 208 nucleons. Parameters are the 
same as for Fig. 1, a small imaginary part e = 0.1 MeV has been added to uj. The 
giant monopole resonance is centred at about 25 MeV because the long-range part 
of the residual interaction has not been taken into account. 

In all the figures showing the multipole strength functions there are two curves. 
The dashed curve corresponds to the strength function associated with the 
response function TZ c ^(tu), while the solid curve corresponds to 1Zl(uo). 

For the monopole response of Fig.3, the two curves are practically coincident. 

For the quadrupole response of Fig.4, the two curves are barely distinguishable, 
while for the octupole response of Figs. 5 and 6, the two curves differ in the 
region of the low-energy octupole resonance, while they are very similar for 
the high-energy component of the giant octupole resonance. 
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Fig. 4. Contribution of the modes (N = ±2, n = 0, ) to the quadrupole strength 
function of the same system of Fig.3. 

The pronounced minimum of the dispersion curve in the top-right panel of 
Fig. 1 at cj ?s 7.7 MeV is at the origin of the increase of the strength function 
given by the solid curve in Figs. 5 and 6. This increase is a collective effect 
due to the pairing interaction. 
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Fig. 5. Contribution of the modes (N = ±1,±3, n = 0) to the octupole strength 
function. The two peaks at oj « 8 and to « 25 MeV correspond to the low and high- 
-energy components of the giant octupole resonance. The energy of the two peaks is 
larger than experiment because the long-range part of the residual interaction has 
not been included in these calculations. The low-energy peak is due to the modes 
N = ±1, while the high-energy peak to the modes N = ±3. 

10 Summary and conclusions 

We have studied the solutions of an improved kinetic equation that includes 
dynamic pairing correlations, as well as finite-size effects, and compared it 
with a more simplified approach based on the constant-A approximation. The 
particle-number symmetry, which is broken in the constant-A approximation, 
is restored by the fluctuations of the imaginary pairing field and all spurious 
contributions are canceled, leading to a correct value of the energy-weighted 
sum rule. The restoration of particle-number conservation is not the only ef- 
fect of the pairing-field fluctuatons since they introduce also the possibility 
of new collective modes of the system, generated by the pairing interaction. 
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Fig. 6. Same as Fig.5, in the region of the low-energy octupole resonance. 



The frequencies of these collective modes are determined by the solution of 
appropriate dispersion relations which are a generalization of the dispersion 
relation of the Bogoliubov-Anderson mode derived by the authors of Ref.[10] 
for a uniform system. In the present paper these modes have been studied 
in detail for a small system like a nucleus, where the pairing gap is small 
compared with the collective frequency. The nuclear mean field has been ap- 
proximated with a spherical square-well potential and it has been found that 
the isoscalar octupole strength function shows an increase in the region of the 
low-energy octupole resonance, which is due to collective pairing effects. 

We conclude that the semiclassical treatment of pairing introduced in [14] and 
improved in the present paper, gives results in good agreement with previous 
quantum approaches and that the self-consistent pairing-field fluctuations con- 
sidered here have two main effects: 

• the fluctuations of the imaginary part are essential for restoring the particle- 
number symmetry which is broken by the constant-A approximation (this 
is already well known from studies of infinite systems); 

• there is a mild collective effect in the octupole response function due to 
the pairing interaction and no similar effect is found in the monopole and 
quadrupole response function. 
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Appendix 



Supplementary condition 

In this Appendix it is shown that only one branch of the solutions of the 
equations of motion (42-45) approximately satisfy the supplementary condi- 
tion (56) required by the Pauli principle for the modes n^O. For this purpose 
we consider the homogeneous system of equations associated with Eqs. (42-45) 
in constant-A approximation: 



u8p™{u) = ujp°*{u) - 2iA<J<(o;) , (165) 

u5p° n d (u) = u n 5p e n v (u), (166) 

ujSk^u) = cvK» - 2A^(w) , (167) 

-iu5K r n (cu) = u^SkKu) . (168) 

The following relation can be easily obtained from these equations, either by 
using all four of them or only the first three: 

(u 2 -uj\- 4A 2 )^» = -2AcVO) • (169) 
It is important to note that, if we put u = w n in this relation, it gives 

= , (170) 



in agreement with the Pauli constraint (56). 

The solutions of the homogeneous system (165-168) are non vanishing only at 
the eigenfrequencies w+(n, I) and w_(n, I), but we have seen that the frequen- 
cies a; + (n, I) are well approximated by ui n and the fact that the condition (56) 
is satisfied when u = Cu n , means that the branch of solutions corresponding to 
w+(n, I) approximately satisfies the supplementary condition (56): 

<KK) « -|«(^) , (171) 

while the same is not true for the solutions corresponding to the frequencies 
w_(n, I). Thus, following Anderson [21] and Valatin [22], we can say that only 
the solutions corresponding to the branch cu+(n, I) ps w„ are acceptable. An 
important point is that the constraint (56) need not be satisfied for all values 
of u, but only at the eigenfrequencies of the system. 
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